!-------------------------------------------------!
!               MODULE GCONTACTS		  !
!-------------------------------------------------!
! This module calculates the green                !
! function of the nano-contacts 		  !
!-------------------------------------------------!

  module Gcontacts
  
    use constants
    use control 

    implicit none
    
    private 
    
!-------------------------------------------------!
! Variables declaration	          		  !					
!-------------------------------------------------!

!-------------------------------------------------!
! Public subroutines	 	        	  !					
!-------------------------------------------------!
    
    public :: G0_contact, G0_wire

  contains

!-------------------------------------------------!
! 	G GENERAL		 		  !					
!-------------------------------------------------!
 
    complex(dp) function  G0_contact(z,epsilon_contact)
    
      implicit none
      complex(dp) :: z
      real(dp) :: epsilon_contact
      
  !    if(contact_shape.EQ.'wire') G0_contact=G0_wire(z,epsilon_contact)

   
    if(WideBandApproximation_logic.EQ.'S')then 
      if(aimag(z).GT.0.0_dp) G0_contact= -cmplx(0.0_dp,1.0_dp)/abs(beta)   !Wide band approximation
      if(aimag(z).LT.0.0_dp) G0_contact= cmplx(0.0_dp,1.0_dp)/abs(beta)   
    else 
!      if(DODContactsFromFile_logic.EQ.'S')then 
!        if(aimag(z).GT.0.0_dp) G0_contact= -Pi*
!        if(aimag(z).LT.0.0_dp) G0_contact= cmplx(0.0_dp,1.0_dp)/abs(beta)   
!      else 
    
    
        G0_contact=G0_wire(z,epsilon_contact)
  !      write(*,*)G0_contact
    endif
    
       
    end function G0_contact
 
!-------------------------------------------------!
! 	G OF THE MONOATOMIC WIRE 		  !					
!-------------------------------------------------!
  
    complex(dp) function  G0_wire(z,epsilon_contact)
    
      implicit none 

      complex(dp) :: z, discriminante
      real(dp) ::  modulo, arg, arg_sqrt, im_z, epsilon_contact, i, j
            
      im_z=aimag(z) 
      
      if(im_z.LT.0.0_dp) z=real(z)-cmplx(0.0_dp,1.0_dp)*im_z   !For the negative part
     
      discriminante=(z-epsilon_contact)**2 - 4.0_dp*beta**2
      modulo=sqrt(Real(discriminante,dp)**2+aimag(discriminante)**2)
      arg=(atan2(aimag(discriminante),Real(discriminante,dp)))
 
      if(aimag(discriminante).LT.0.0_dp) arg=arg + 2.0_dp*pi
      
      arg_sqrt=arg/2.0_dp       
           	
      G0_wire=(z-epsilon_contact) - sqrt(modulo)*exp(cmplx(0.0,1.0,dp)*arg_sqrt)
      G0_wire=G0_wire/(2.0_dp*beta**2)

      if(im_z.LT.0.0_dp) G0_wire=real(G0_wire)-cmplx(0.0_dp,1.0_dp)*aimag(G0_wire)   !For the negative part
      z=real(z)+cmplx(0.0_dp,1.0_dp)*im_z
               
	   
    end function G0_wire

!-------------------------------------------------!
! Here you can calculate others G functions       !					
!-------------------------------------------------!

  end	
